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We study the Effective Field Theory of Large Scale Structure for cosmic density and momentum 
fields. We show that the finite part of the two-loop calculation and its counterterms introduce an 
apparent scale dependence for the leading order parameter cj of the EFT starting at fc = 0.1 /iMpc“^. 

These terms limit the range over which one can trust the one-loop EFT calculation at the 1% level 
to fe < 0.1 /iMpc“^ at redshift z = Q. We construct a well motivated one parameter ansatz to fix 
the relative size of the one- and two-loop counterterms using their high-fc sensitivity. Although this 
one parameter model is a very restrictive choice for the counterterms, it explains the apparent scale 
dependence of Cg seen in simulations. It is also able to capture the scale dependence of the density 
power spectrum up to fc « 0.3 hMpc“^ at the 1% level at redshift 2 = 0. Considering a simple 
scheme for the resummation of large scale motions, we find that the two loop calculation reduces the 
need for this IR-resummation at fc < 0.2 /iMpc“^. Finally, we extend our calculation to momentum 
statistics and show that the same one parameter model can also describe density-momentum and 
momentum-momentum statistics. 


I. INTRODUCTION 

The development of an Effective Theory for Large Scale Structure (EFT of LSS) [U [5] has lead to a resurgence 
of interest in perturbative approaches to study the development of structure in our Universe. Although Standard 
Perturbation Theory (SPT, see e.g. [S]) has allowed invaluable insights in the physics of LSS, it was soon realized 
that in order to extend the validity of the theory down to smaller scales, one needs do go beyond SPT. Various 
approaches that can be found in the literature, e.g. Hd, focus on the resummation of higher order contributions in 
order to achieve an accurate description of non-linear data up to large wavenumbers. Perturbation theory will never 
be able to capture the small scale dynamics even after complete resummation (as an illustration see 0)- This fact 
limits the applicability of resummation results. The EFT of LSS aims at extending SPT through the modeling the 
effects of small scale dynamics on larger scales. Based on the same principles that high-energy physics community has 
been exploiting for decades, the EFT framework allows to describe perturbatively the evolution of long wavelength 
modes while systematically taking into account the impact that short wavelength modes can have. The power of 
this approach lies in the fact that it is not necessary to explicitly solve the non-linear small-scale dynamics, which, 
however, comes at the cost of introducing parameters that are not determined by the theory itself. Furthermore, the 
EFT approach allows to overcome some conceptual shortcomings of SPT. 

Since the original papers mm. many aspects of the EFT of LSS have been explored in the literature. At the 
one-loop level [SHU] have made progress, while in min] attempts were made to tackle the EFT of LSS at the 
two-loop order. These references mainly focussed on the two-point functions, while usmsi considered the matter 
bispectrum and EOllll] looked at non-Gaussian effects. The Lagrangian space formulation of the EFT of LSS has 
been studied in [22[ |23] and finally, aspects of bias and baryonic effects have been considered in [201 1241130] . 

Despite the power of the EFT approach, we have to deal with the presence of free parameters. For the power 
spectrum at the one-loop level, one such parameter is introduced. For the bispectrum, three additional parameters 
are necessary and for the two-loop power spectrum an even larger number of free parameters would have to be 
consideredjj This, however, means that the theoretical description stops being predictive and three- and four-point 
functions would have to be considered in order to determine the values of all parameters through a comparison with 
simulations or observations. 

The role of the free parameters and the corresponding counterterms is to incorporate the effects of the small scales 
into the theory. Our aim is to achieve exactly this in a systematic way for the power spectrum at the two-loop level. 
The approach that we follow does, however, avoid an unmanageable number of free parameters by making a well 
motivated ansatz. 
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^ Consider e.g. Chiral Perturbation Theory where in the strong sector at the leading order there are two, at the next higher order 12 and 
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In this paper we revisit the two-loop Eulerian power spectrum calculation and compare the results to our own 
set of numerical simulations. We compare results at the level of the power spectra as was done in [Ml Hz]. In two 
companion papers EH [32] we compare perturbation theory with the results of numerical simulations for the same 
initial conditions. This is a more stringent test than what is presented here. Our goal in this paper is to reproduce the 
comparison method used in the literature and try to relate the result to what we see in the more detailed comparison. 
We will find that in terms of the maximum k where the perturbative calculation can be trusted both results agree. 

This paper is organized as follows. After a brief review of the EFT of LSS, we consider in more detail the UV 
sensitivity of the one- and two-loop integrals in Sec. jllBj and jll C| From the UV sensitivity, we derive our ansatz for 
the two-loop counterterms in Sec. II Dj and an even simpler procedure for the counterterms is discussed in Sec. HE In 
Sec, jlllj we compare our approach with numerical simulations and present our results. Also, we discuss the two-point 
correlations functions that involve momentum. 


II. THE EFT OF LSS 


In the EFT of LSS one sets to solve perturbatively the following equations: 

drS + di[{l + i5)w*] = divJ- , 


drv^ + nv^ + = - a-r*- 

ap 


( 1 ) 




.s. 


These equations differ from those of SPT |3] due to the addition of new source term, u* in the continuity equation 
and a stress tensor source r®-! in the Euler equation. These sources arise from small scales, where the perturbative 
solution of SPT is not applicable. In the EFT of LSS they have to be modeled as they arise from modes that are 
outside the range of applicability of the theory and thus result in the introduction of free parameters. The EFT of 
LSS provides an organizing framework for how to model these sources, providing a list of terms with their associated 
free parameters that need to be introduced to achieve a desired accuracy. 

For simplicity, in the discussion that follows we concentrate on the stresses that appear in the Euler equation. In 
this paper we will not consider velocity statistics, but only statistics involving the density and the momentum. In 
such case it suffices to discuss the stresses in the Euler equation as the effects from u® in the statistics we will consider 
can be mimicked by changing r®-^. In any case, all the conceptual points we will make below are applicable to both 
T®^ and u®. 

The T®-^ stresses come in two different forms. Some of these stresses can be computed in terms of the perturbative 
solution, others cannot. For the latter one only has a model for the statistical properties of those stresses. It is 
convenient to decompose the velocity field into its gradient and curl pieces. At the order we will work in this paper 
only the gradient component will be relevant, thus the stresses we need to model only enter through a scalar quantity: 


Te = -di 


1 

ap 


djT^ 


^det 


+ T‘, 


stoch 


The deterministic part of the stresses Tg can be modelled perturbatively. In the EFT we write schematically 

rr = 


( 2 ) 


(3) 


The deterministic part of the stresses is a local function of the perturbative solution, and we have used the equivalence 
principle to demand that it can only depend on second derivatives of the gravitational potential (higher spatial deriva¬ 
tives and time derivatives can also appear). We have introduced 4> = so that didj^ is dimensionless 

and = 6. For the stochastic part, all we can do is model the statistical properties of 

In the EFT of LSS is modeled as a power series in didjcj) and its spatial and time derivatives. In addition 
to the equivalence principle, mass and momentum conservation restrict the form of both and of the statistical 
properties of In particular in Fourier space Tg®*(/c) needs to go to zero at least as faster than the density 

when k ^ 0 and the power spectrum of should go to zero at least as k'^. 

To calculate the one-loop power spectra in ACDM, only the lowest order piece of r^®* is relevant. It is given by 


^det I 


ILO 


= —(i^A(5(i) = —d^AA(/)(i) 


(4) 


where 5(i) is the linear solution of perturbation theory. In this formulation, because r^®* acts as a source in the 
equations of motion, the time dependence of cP will affect the results. In particular it will be relevant to determine 
the relative sizes of the corrections in the different two point functions involving 6 and 9. 
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FIG. 1. Diagrams for the tree level, one- and two-loop expressions of the SPT power spectrum. 


The case of the one loop bispectrum has already been considered in the literature HHiiin]. In that case the second 
order counterterms are needed. This introduces three additional parameters for the spatial structure of One can 
write: 

= -^^^^['^( 1 ) + <^( 2 )] - (5) 

with 

( 6 ) 

In principle, d, ei, 62 and 63 could be fixed by measuring both the power spectra and bispectrum. In practice however, 
with current simulations there are signihcant degeneracies among these different parameters. Making an ansatz for 
the ratios, scaling all counterterms by the same amplitude and fitting for this overall amplitude parameter, seems 
good enough to explain simulation measurements [18] . 

In this paper we are interested in performing a two-loop calculation for the power spectrum and thus we would have 
to model the stresses up to third order in the helds. Modeling these terms will increase the number of parameters even 
further. At the level of the two point function however, some of these parameters will be degenerate. In principle, 
one could disentangle all the new parameters comparing the predictions with the four point function measured from 
simulations. In practice the necessary signal to noise ratio to do this is probably not available in the current generation 
of simulations and a simple ansatz for the ratios of amplitudes of the various terms could be good enough. In any 
case, in this paper we will only compare results against measurements of the two point function and thus we will not 
have enough information to distinguish all the parameters. Furthermore, in this type of exercise one runs the risk 
of overfitting the power spectra simply because one is introducing too many additional free parameters. In order to 
avoid this, one should compare the results of perturbative calculations with simulations at the level of the fields as was 
done in |31j for the Lagrangian displacement and in |32j for the density. In this paper we will adopt a simple ansatz 
for the size of the various counterterms and only keep one overall free amplitude as a parameter. We will discuss this 
in the next sections. 


A. Perturbative solution and counterterms 

In Standard Perturbation Theory (SPT, for a review see[3|) the perturbative solution of the equations of motion 
has the following structure, 

S = (5(1) -I- 5(2) + 5(3) -I- 5(4) -I- 5(5) -!-••• (7) 

where 5(„) depends on the initial conditions to the n-th power and we have only written terms relevant for the two loop 
calculation of the two point function. When computing the power spectrum, one considers the averages of {S(^n)S{m))- 
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FIG. 2. Diagrams that are regularized in our approach. The dashed loops are the ones where the momenta is large and are 
fixed by a counterterm. 

At the tree level, the only possible order is n + m = 2 and therefore n = m = 1. For one loop, n + m = 4 so the two 
possible terms are the mixed term between i5(3) and i5(i) or the square of 5(2p At two loops we have n + m = 6 so the 
options are 1-5, 2-4 or 3-3. In SPT one writes the perturbative solutions a^ 

<5(n)(fc) = [ ■■■ [ i2TrfS‘^^\qi + ...qn-k)Fn{qi,...,qn)So{qi)...6o{qn), (8) 

•Jqi J Qn 

where Jq stands for the initial density fluctuations. The different contributions to the power spectrum computation 
can be represented using the diagrams in Fig. and combine to the power spectrum as 

P 55 = Pll + 2Fi3 + P 22 + 2Fi5 + 2 P 24 + T 33 -I + ^33-11 • (9) 

The explicit expressions for the constituent power spectra are given in App. The integrals for the one- and two- 
loop contributions to the above expression bear some UV-sensitivity or can be even divergent for certain input power 
spectra. The EFT provides a framework in which these UV-sensitivities can be addressed and regularized with the 
corresponding counterterms. That is to say, that the EFT counterterms provided by the stress tensor and its time 
dependence in Eq. Q should be able to capture and correct the UV-sensitivity of the SPT expression. 

The equations of motion Q only have quadratic non-linearities, so vertices in diagrams should only be cubic. That 
is to say, the Fn kernels in the diagrams we showed in Fig. [l] are effective time integrated diagrams that can be 
constructed by having multiple cubic vertices joined by propagators (or Green’s functions) [1]. In the EFT, there are 
additional diagrams due to the introduction of counterterms, or sources in the equations of motion. 

The first thing one notices is that some of the two-loop diagrams contain inside of them a subdiagram that looks 
like a one-loop diagram. The EFT procedure amounts to adding a counterterm that corrects the mistakes introduced 
when a high momenta is running in a loop. This is schematically shown in Figure If at least some of the two- 
loop diagrams contain pieces that look like one-loop ones, then the same counterterms that have fixed the one-loop 
subdiagram would fix the two-loop ones. To accomplish this, one would need to solve the equations of motions with 
the one-loop counterterm as a source to obtain a solution linear in the amplitude of the counterterm but up to cubic 
in the initial conditions 5q. Because the counterterm will be a source in the equations of motion acting over time, 
carrying out this calculation would require specifying the time dependence of the one-loop counterterm. This program 
was carried out in m as well as in the case of the bispectrum [181 US] ■ 

But even for the diagrams that naively look like those in the one-loop calculation, putting the one-loop counterterm 
into the equations of motion does not necessarily fix all the loops correctly. The point is simple: in the one-loop case, 
the diagram is computing the effect of a short mode that evolves in a linear long wavelength background. Thus, the 
time evolution of this background is given by the linear growth factor. In some of the two-loop diagrams the short 
modes in the loop are evolving in a background that is quadratic or cubic and thus the details of their evolution 
and the value of the counterterm could not be the same. This fact was already noted in the one-loop bispectrum 


Note our shorthand notation for the integral and measure 


Ir-I 


d^q 


(27r)3 ■ 

Furthermore we will express momenta with respect to the external momentum as q 2 = r 2 k and qi = rik. 













5 



FIG. 3. Effect of changing the cut-off from Ah = 5 hMpc“^ to Ai = 1 hMpc“^ for the one and two loop calculations normalized 
by k^P. Left panel: Contributions from the low-high and high-high terms (single- and double-hard). The mixed term clearly 
dominates the k^P part and also the deviations from this behavior. Right panel: Contributions from the separate diagrams. 
At the one loop level P 13 leads to a k^P contributions, whereas the contribution from P 22 is suppressed. P 15 dominates the 
k^P part but for the deviations from this scaling, there is a cancellation between P 15 , P 33-11 and P 24 . Like P 22 in the one loop 
case, the term arising from P 33-1 is suppressed. 


calculation, where it was shown that the counterterms coming from the time evolution of the linear counterterm are 
not able to capture the UV-sensitivity of the SPT loops [18]. 

Of course, in addition to the terms that derive from the linear counterterm through the equations of motion, there 
are those that arise from the new quadratic and cubic contributions to the stress tensor. Once all of these counterterms 
are included, one has sufficient freedom to correct all UV mistakes at this order. The entire set of counterterms could 
be fixed by studying the power spectrum, bispectrum and trispectrum. 

In [T7], the first two-loop calculation in the EFTofLSS, only the counterterms that follow from the leading order one 
were kept. Thus, the one- and two-loop counterterms depended on only one free parameter (and its time dependence). 
This was done mainly for simplicity, as one could not fit multiple parameters simultaneously from the available power 
spectrum data. Here, we will take a similar strategy, in that we will also study a one-parameter family of counterterms, 
but we will fix them in a different way. Because in both, our calculation and in El. one is using an ansatz for the 
two-loop counterterm, one should recognize that in all generality its amplitude could be somewhat different than the 
one being calculated. 


B. UV-sensitivity at one loop 


The counterterms introduced in the EFT are there to model the effects of the small scale dynamics on larger scale 
modes. Thus, a place to look for an ansatz for the relative sizes of the EFT counterterms could be to study the effect 
of a shell of power at high loop momenta computed in SPT. We can fix the ratio between the various counterterms in 
the EFT to the one given by this ansatz but leave the amplitude of the small scale power in the shell as one overall 
free parameter. This ansatz makes the final results insensitive to the small scale power in the SPT calculation. 

For this purpose, we start by computing the contribution of a shell in momentum space between Ah = 5 /iMpc”^ 
and Ai = 1 /iMpc”'^ to the one-loop power spectrum in SPT. This choice is somewhat arbitrary, but provides us with 
a sufficiently significant change to see the effects and furthermore the lower limit is sufficiently far away from the 
scales of interest fc « 0.1 hMpc”^ to warrant a separation of scales. We call this contribution to the one-loop power 
spectrum = 2Pi^ + P 22 and the results are shown in Fig. We recover the standard result, that for ACDM 

at small wavenumbers the P 43 contribution dominates and scales as k^Pn- The subdominant PI 2 contribution scales 
as In this language, what is usually called the correction in the EFT is nothing other than the functional form 
of the fc —>■ 0 limit of P^/k’^Pw. 

The value of P^oop depends not only on the amplitude of the power added on the shell but also on the position 
of the shell. In our ansatz for the counterterms, we can also use analytic expressions obtained in the limit that the 
momentum of the shell (qi) is much bigger than the momentum of interest k qi. We will call this limit Pil^^, 
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basically moving the shell to infinity. For the leading and sub-leading contributions we obtain: 


piji—foo _I,2 p 

13 - ^ ^ 630 Jg 


61 f Pi<l) , jAp 2 f P{q) 


/2 = -^n -2 


( 10 ) 


The leading contribution is proportional to the high-g contribution to the one dimensional displacement dispersion 
f dq P(g) and the shell power will thus be = —61/210A:^(t^ sh-P’n- shell under consideration 

here, we have gj, = 0.68 h~^Mpc^. It is interesting to note that the square of this coefficient has a factor 100 
stronger cut-off dependence than the coefficient of the subleading kf^P contribution, which makes sense since the 
integral of the latter is suppressed by two additional powers of q in the UV. Our strategy will be to add the effect of 
this shell of power computed up to two-loops to the standard SPT results with a free parameter. 

This approach is equivalent to the EFT, where the effect of a high-A: shell can be captured by the leading counterterm 
Cg (it should be noted that Cg is not equal to the parameter tP in Eq. as Cg is the result of a time integral over the 
Green’s function and (P). In terms of the standard notation in the literature, 

Pctr.iioop = -2k‘^clPii (11) 

and at the level of the density field it corresponds to the term 5(i) = —P^k'^Spy It is common practice to fix the 
coefficient Cg after the full one-loop calculation has been subtracted from the data and we will adopt this convention 
here. Thus, the true coefficient K 2 of the k'^Pn part of the low-fc limit of the data is fixed to be 

-2c^ (12) 

Consequently, the number Cg effectively contains all the higher order loop contributions to K 2 , their counterterms and 
the true small scale contribution. In particular, no higher loop contributions to K 2 should be calculated and to the 
extend that such terms are present in higher order calculations, they should be removed. 

Numerically, we will find below in accordance with previous studies that Cg is a positive number of order 1 h~^MpP. 
This means that we are effectively increasing the power in a high-fc shell in perturbation theory, but the effect is an 
enhanced large scale damping of the non-linear power spectrum. 


C. UV-sensitivity at two loops 

We now evaluat^ the total two-loop power spectrum as well as its constituent pieces for two different cut-offs 
Ai = 1 /iMpc“3 and Ah = 5 /iMpc”^. As we show in Fig. adding this shell of power primarily affects the k^Pn 
coefficient. This piece should be absorbed by the counterterm that was already introduced at one loop in the previous 
section. The only relevant parts are the deviations from the k^Pn behavior for fc > 0.1 /iMpc”^. These are the 
deviations that we want to capture. These non-trivial pieces should be captured by the two-loop counterterms. 

In contrast to the one-loop calculation, we now have two momenta that are integrated over and thus we have to 
distinguish two cases: i) both loop momenta are large (both loop momenta in the in the high-A; shell, double-hard) 
or ii) only one loop momentum is large with respect to the other momenta in the problem (one momentum in the 
high-A: shell, single-hard). The left panel of Fig. shows these two contributions separately. We immediately see 
that the double-hard limit is basically degenerate with the k^Pn behavior for all the ks of interest and thus is not 
very relevant for our calculation. There is a slight upturn for high wavenumbers that we will discuss in more detail 
below. The single-hard contribution also has a k^Pn part, in which we are not interested, but beyond this it has the 
interesting new scale dependence that should be captured by new counterterms. This motivates us to consider the 
single-hard limits of the two loop calculation. 

The terms leading to the k^Pn contribution in the shell calculation are also present in the finite part of the two-loop 
calculation, actually governing its low-A; behaviour. If this contribution was kept in the final calculation, it would 
change the value of the parameter of the one-loop counterterm Cg in Eq. ([n ) that was introduced to regularize the ^3 
contribution. As we stated in the previous section, our strategy is to fix this number after the one-loop calculation. 
We thus decide to remove the k^Pn term from the finite part of the two-loop calculation. This can be done by either 


^ The numerical integrals for the two-loop expressions are performed with the CUBA libraries m SUAVE routine employing the IR-safe 
integrand [3II6|. 
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calculating the limit of P 15 as we do in Eq. (C 8 ) or by fitting the very low-fc limit of the numerical calculation. We 
will denote the two-loop calculation, from which the degenerate part has been removed as P 2 ioop- More generally, all 
terms that have been corrected for degeneracies with lower order counterterms will be decorated with an overbarj^ 
Let us now come back to the limits of the two-loop calculation. Before we discuss the single-hard limits that we 
deemed responsible for the new counterterms at two loops, we discuss one double hard limit, namely the one of P 33 _//, 
that leads to 


nil,<12^00 

^ 33-11 



(13) 


The counterterm for this term is automatically included once the square of the leading order term at the field level 
(^( 1 )^( 1 )) (X is considered. There is also a piece of P 15 that has this same structure. 

As before, we can obtain analytical formulas when the shell is taken to be at infinite momenta. Because we only 
care about the piece that does not look like fc^Pn, we focus on the case when one of the two momenta is running in 
the loops is in the shell {qi) while the other one remains finite ( 92 ), *-e., the single-hard limit. Fig. [^shows that the 
total effect from the two loop terms arises from a cancellation between upturns in P 15 and P 33 -J 1 and a downturn in 
P24: while P 33-1 is basically flat. The overall effect is a residual suppression of power on small scales. We would not 
have needed to consider an explicit shell, but could have considered all the limits of the two loop calculation (as we 
do in Appendix]^ to see which terms will lead to relevant counterterms. We considered the approach presented here 
more pedagogical. 

Leaving the double-hard limits and the suppressed single-hard limits for discussion in Appendix let us give here 
only the relevant terms. The formulae for these single-hard limits read: 


pqi—>-oo _ 2 


X 5 X 3y ;^271350^08000r5 ] “2^2(2266005 - 33470730r^ -b 187902172r^ - 9879110r^ -b 1167375r^) 


-bl5(-l -b r|)3(-151067 - 451074ri -b 77825r|) log 


1 -b r2 
ll-rzl 


>PiQ 2 )crdk P{k) , 


(14) 


= 4 X 3 X 3fcV^ f 
J q 2 

xP2(g2,fe - g2)P(|fe - g2|)P(<?2) , 
P!3-n = -2aj-^k^P33{k) , 


(-32879/42 + r|(6176 - 48096/4^) -b 32r|/i2(1117 -b 1503/tl) -b r2(-25933 -b 16892/4^)) 


4074840r2(l-b rl - 2r2/42) 


(15) 

(16) 


where we defined = qi/k and k ■ qi = ^ikqi. Note that all three single hard limits are proportional to the high-g 
contribution to i.e., the small scale displacement dispersion. As for P33_!1}Y'°° the counterterm for P^l^^ is 

automatically included once the leading order counterterm at the held level is correlated with the third order held, 
leading to oc — 2 clk'^Pi 3 . The term corresponds to the UV-sensitivity of the bispectrum term B 114 

calculated in [18] and is thus hxed by the corresponding counterterms. The term corresponds to a contraction 

of the UV-limit in the trispectrum term Tins and should thus be regularized by the corresponding counterterm. 

Slightly problematically, the part of the integral is log-sensitive in 52 for q 2 ^ k 


= -2 X 5 X 3 X 3ajk^P{k) y yg|lp(g 2 ) . (17) 

This limit is proportional to k^Pn and thus completely degenerate with the leading order counterterm c^. It is this 
limit that leads to the offset in the single-hard limit in the left panel of Fig. It is continuous with the limit in 
which in a first step qi and 52 become large (the double hard limit of Eq. |C8| ) and then one momentum is taken to 
be smaller than the other one 52 ^ 9i- To the extend that we are not interested in this contribution we define a new 
limit from which this term is removed at the integrand level 


pqpoo ^ + 2 X 5 X 3 X 3ajk^P{k) 


120424 

45147375 


^( 92 ) . 


( 18 ) 


^ In 1171 the subtraction of the two-loop contribution proportional to k^Pii was accomplished by introducing the parameter P 











D. Ansatz for the two-loop counterterm 


Our ansatz to have a one parameter counterterm is to use: 

P^‘' = a«op + Atoop), (19) 

where a is an overall free parameter that can be set by looking at the piece of the power spectrum that scales as 
k'^Pii at very low k and -P|ioop is the two-loop power spectrum obtained from the shell calculation with the degenerate 
k^Pii piece removed. As we have seen above, all the relevant limits at one- and two-loop are proportional to cr^ 
that we can relate our ansatz to the usual EFT language by setting a = 210Cg/(61cr^). 

More explicitly, for the counterterm at the two loop level we will consider 


such 


Petr ^ 

= P 








^33-11 


ctr,lloop 


P 


ctr,21oop 


( 20 ) 

All the terms in the right hand side of Eq. (20 1 are proportional to and in fact these are the only terms proportional 
to tr^. Hence, our ansatz is nothing other than choosing the value of tr^ by matching the low k behavior of the power 
spectrum from simulations to the k^Pn template. Since this is an important point, let us repeat again the basic idea 
of our approach. The relation between cl and cr^ that is found at the one-loop level in Eq. (12) is used in order to 


cure the UV sensitivity of the two-loop integrals. Effectively, for all occurrences of the problematic tr^ in the two-loop 
integrals we add a cl counterterm. This is what is shown in Eq. (20) and we end up with a one-parameter model for 
the UV sensitive parts of the one- and two-loop integrals. 

Finally, note that the standard IR cancellation when q 2 k qi still happens among the qi-limits computed 
above: P 15 , P 24 and P^^-jj. In this case the l/g| motion contributions cancel and only long wavelength tides survive: 


pgi->-00,(32- 

^ 33-11 


= Pik) [ f 

Jqi -'q2 


3538 


-b 


61 1 


2p|1^OO,92^0 ^ f f 

J qi J q2 

2pq^l^oo,,32^0 ^ 


99225r2 1890 


1361863 61 1 


qi •Jq2 


5942475r2 945 


902354 ^ 61 1 

'4729725r2 1890 ^ 


^0 _|_ >- 00 , 92—^0 >-0 _ P(k) 


P(<7i)P(<72) , 

P{qi)Piq2) 

Piqi)P{q2) 

12670991 


qi 'J q2 


27810783r^ 


Pi.qi)Pi.q2) ■ 


( 21 ) 

( 22 ) 

(23) 

(24) 


E. An even simpler ansatz 


Finally, we could consider what is perhaps the simplest ansatz of all. Just as in the case of the stresses parametrized 
by Ts in equation ([^, one can parametrize the counterterms in such a way that one of the terms is just proportional 
to the density computed in SPT. That is the terms relevant for the two-loop calculation at the level of the fields could 
be written as: 

5 = J( 1 ) -b J(2) -b i 5(3) -b (5(4) -b (5(5) — /^A((5(i) -b (5(2) -b (5(3)) -b ... , (25) 

where the ellipsis account for the terms arising from other quadratic and cubic counterterms. An extremely simple 
ansatz is to set those additional terms to zero. This would lead to the following expression: 


P 


ctr,simple 


= -kH 


2;2 


Pi + Plloo 


(26) 


As we discussed earlier, neither of these ansatzes is expected to be perfect, and nothing short of fixing all the 
counterterms by studying the three- and four-point functions or projections at the field level would be perfect. The 
philosophy of this paper is to write down examples which are expected to have roughly the right size and use those 
to asses how big these terms are expected to be while keeping in mind the uncertainty in their size. Fig. [^compares 
the two formulae proposed in this section. They are in reasonable agreement in terms of both the expected shape and 
the size of the correction. One could definitely argue that Psimpie is perhaps too simplified as one is ignoring effects 
that we know are there and are furthermore comparable to those being included. Our P^tr defined in Eq. (20) has all 
the relevant terms included, although perhaps some of their relative amplitudes are not correct in detail. 
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FIG. 4. Comparison between the two-loop connterterm deduced from the divergencies and the one-loop power spectrum 
weighted by wavenumber squared. We see that the explicit calculation of the two loop counterterms Pctr, 2 ioop is proportional 
to the naive estimate k^Puoop- 


F. Relative size of the corrections 

Before comparing with simulations, we can take our model for the power spectrum and calculate the sizes of the 
different terms. This comparison will allow us to estimate how well we expect our formulae to agree with simulations 
and estimate reach of perturbation theory. In particular, now that we have an estimate of the two-loop terms and 
their associated counterterms we can ask when they make a difference relative to the one-loop terms and ask over 
what range of ks it would be safe to fit for when doing a one-loop calculation only. Our full two-loop EFT power 
spectrum is 


P = Pii + Plloop P .fctr.lloop .f^21oop P .fctr,21oop- 


(27) 


We can now compute two quantities: 


- 1 = 


P 

P — Pll — Plloop 

2fc2Pii 


.f^lloop P .fctr.lloop P P2\oop T .Pctr,21oop 


= c!- 


P 


21oop 


-Pll 

■ -Pctr,21oop 


2fc2Pi 


(28) 


The first of these two quantities indicates the size of the various terms as contributions to the power spectrum, the 
second indicates the relative correction they would make to a fit of after subtracting the explicit one loop SPT 
calculation from the data. We show these quantities in Figure]^ The left panel shows that both Aioop and 2ioop 
make roughly a 5% correction to the power around k = 0.2 /iMpc”^. Given that Pctr,2ioop is uncertain because 
we have not used three- and four point function measurements to obtain its amplitude but only have an ansatz, it 
is difficult to imagine that one could be more accurate than about one percent on these scales. The counterterm 
is relatively steep, so even though it contributes 5% around k = 0.2 hMpc”^ at fc = 0.5 /iMpc”^ is makes an 
order unity contribution. In the same panel, we also show the effect of the c^k^Pn correction, which is at the sub 
percent level for the wavenumbers considered here. Note however, that the coefficient of this term should be fitted 
independently, since it has to capture the subleading UV sensitivities in P13 and P15 that we have neglected so far. 
We also estimate the three loop counterterm at the basis of our most simple counterterm ansatz, i.e., we consider it 
to be given by —k‘^P 2 ioop- This term leads to percent level corrections at k = 0.3 hMpc”^, so we should be worried 
about similarly large corrections from the three loop calculation for even larger scales. Finally, we overplot the size 
of the stochastic term estimated in |32]. Given that it leads to percent corrections at k = 0.25 /iMpc”^, we should 
not expect any perturbative approach to match the full power spectrum to a better accuracy than this. Actually, the 
perturbative/deterministic calculation performed here should describe the non-linear power from which the stochastic 
part has been removed. 

As we show in the right panel of Fig. when fitting for c^, the combination Aioop + -Pctr,2ioop changes by about 
50% between k = 0.05 hMpc”^ and k = 0.20 hMpc”^. About half of this change is from the finite part of the two-loop 
calculation and half from the counterterms. These two corrections are of the same amplitude at fc = 0.18 /iMpc”^. 
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FIG. 5. Left panel: Estimate of the size of the corrections arising from various contributions to the two-loop calculation. 
The finite part of the two-loop calculation leads to percent level corrections at fc = 0.1 /iMpc“^. We also show the corrections 
from the square of the speed of sound term which is suppressed over the range considered here. The size of the three- 

loop counterterm can be estimated as 0{1) x k^P 2 \oop and leads to percent level contributions at k — 0.3 hlVIpc”^. We also 
show the estimate for the stochastic part of the total power spectrum from [32] which leads to percent level corrections at 
k = 0.25 hMpc”'^. Right panel: Estimator for the leading EFT coefficient cf. The model is evaluated for = 0.98 h“^Mpc^ 
and the gray band shows the effect of a 10% change in this value. Note that at fc = 0.2 ZiMpc”^ the one-loop counterterm 
and the two-loop correction are of the same order. The two-loop term leads to a considerable scale dependence of for 
k > 0.07 /iMpc“'^. 


Due to the presence of these corrections, a measurement of Cg without consideration of the two-loop terms is not 
possible for wavenumbers exceeding k = 0.07 /iMpc”^. Besides the broadband upturn, there are also considerable 
wiggles from the Baryon Acoustic Oscillations (BAO) in the finite part of the two loop calculation. 

At this point it is perhaps instructive to write an equation relating the change in the inferred value of (Ac^) to 
changes or errors in the power spectrum (AP): 


Ac 


2 

S 


AP 1 


AP/P / k 
2% U-1 hMpc"^ 


h-^Mpc^ 


AP/P / k 
0.2% V0.03 /iMpc"^ 


h-^Mpc^ . 


(29) 


For values of around 1 h“^Mpc^ and a measurement at fc ~ 0.1 /iMpc”^, an accurate measurement of requires 
one to model all other contributions to the power spectrum at the sub-percent level. Besides that, the statistical error 
should also be at this level. State of the art simulation codes and reasonable simulation volumes can deliver this level 
of accuracy and precision. However, our estimates above show that at this scale, one needs to include the two-loop 
terms. If one goes to fc ^ 0.03 /iMpc”^ higher loop contributions are negligible, but as we discuss later the required 
10 “^ level numerical precision might be challenging. 


III. COMPARISON WITH SIMULATIONS 

As a benchmark for the performance of the perturbation theory we employ a suite of dark matter only simulations 
of the WMAP7 cosmology [M] (17^ = 0.272, Da = 0.728, rig = 0.967, erg = 0.81). We have run 16 simulations with 
a box length of 1500 /i“^Mpc (L simulation) and also one realization of a smaller size, higher resolution box with 
500 /i“^Mpc box length (M simulation). The simulations are initialized with the second order Lagrangian Perturbation 
Theory code 2LPT |35] at redshift Zi = 99 and the 1024^ particles are subsequently evolved using GADGET2 |36] to 
redshift z = 0. For more details on the simulations and some convergence tests see m- 
For the speed of sound in the one-loop EFT calculation we employ the following estimator 


^2 7^11 7\ioop 

“ WlAi 


(30) 


where Pni is the power spectrum from the simulations. In Fig. we show the measurements at redshifts z = 0, 0.5,1, 2 
from our simulations. The data clearly show a scale dependence with significant deviations from the low-fc limit at 
higher wavenumbers. There are also distinct BAO wiggles in the measurement that have been noted in the literature 
m- We have corrected the data for 2 x 10 level deviations in the linear growth factor, that are likely connected to 
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FIG. 6. One loop Cg estimator at 2 = 0,0.5,1,2 from top left to bottom right before IR resummation. At all redshifts, we 
see clear evidence for a running of Cg that is described by the scale dependence of the two loop correction and the two loop 
counterterm. 


the integration accuracy in GADGET and would lead to a low-fc upturn in this figure. Furthermore, we have cancelled 
the leading order cosmic variance, by actually considering the ratio of non-linear power spectrum and linear (initial) 
power spectrum measured in the simulations. The one-loop EFT model (horizontal magenta dashed line) fails to 
describe the data for k > 0.07 ZiMpc”^, but the two-loop corrections can explain the residual scale dependence. We 
find Cg « 0.98 h~^Mpc^ and interpret the difference from previous measurements Cg « 1.6 h“^Mpc^ m extracted 
from the k = 0.15 — 0.25 /iMpc”^ range as resulting from the two loop contributions. It is also worth noting that the 
two loop calculation is already doing a very good job at tracking the BAO oscillations, at least for k < 0.2 /iMpc”^. 
The calculation based on the UV-limits assumes a time dependence of Cg that matches the one of the SPT term 
that it is regularizing, i.e., D^{a). We are using this time dependence to scale our z = 0 fit to higher redshifts 
and find very good performance both for the small wavenumber behaviour as well as the scale dependence at higher 
wavenumbers. All of the redshifts show slightly low datapoints at fc = 0.03 /iMpc”^ and k = 0.045 hMpc”^, that 
spoil a nice asymptotic behaviour at low wavenumbers that one would expect in the EFT. As we describe in App. 
this systematic effect goes away if the theory is calculated on the simulation grid, effectively using the same modes 
that are present in the simulations. 

The discussion of the relative difference between simulation and analytic calculation for the power spectrum itself 
will be defe rred until we discuss the IR-resummation below in Sec. IIIB[ but the anxious reader might want to look 
at Fig. 10 The two-loop calculation agrees with the data at the sub-percent level all the way to fc = 0.3 /iMpc”^ at 
z = 0. 


A. Time Derivative and momentum correlators at two loops 


So far we have concentrated at a single observable, the density power spectrum, we would now like to extend 
the calculation to momentum statistics. This extension is motivated by the fact that the momentum statistics are 
sensitive to the time dependence of the speed of sound and more sensitive to loop corrections than the density power 
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FIG. 7. Constraints on the speed of sound and its time dependence from the momentum-density cross power spectrum. 


spectrum itself. We consider momentum /r = —V • [(1 -I- instead of velocity, since the latter is only defined at the 
particle locations in the simulation and thus the only quantity that can be reliably measured. To compute correlations 
involving momenta we then use the continuity equation: 

fi{k,a) = S'{k,a) . (31) 

For the SPT predictions of the momentum-density cross power and momentum-momentum auto power we have 


=fn Pii + 2 (2 Pi3 + P 22 ) + 3 (2Pi5 + 2 P 24 + P 33 ) , 




Pii + 2 (3Pi3 + 2 P 22 ) + IOP 15 + I 6 P 24 + 9 P 33 


(32) 

(33) 


Note that while P^^ is IR-safe, i.e., the cancellation of IR modes between P 13 and P 22 still happens (and similarly for 
the two loop contribution), this is not the case for P^fj,. The density-momentum cross correlation is nothing than the 
time derivative of the density power spectrum 

(<5(fc),Mfc')) = {m,S\k')) = {27Tfd^^\k + k')^P^^s{k) . (34) 

Thus, if we manage to capture the time dependence of the perturbative corrections to the density power spectrum, we 
should also be able to describe the density-momentum cross power spectrum at intermediate times. We immediately 
see that the importance of loop corrections is enhanced with respect to the density-density power spectrum and thus 
expect a steeper running of the speed of sound corrections. 

Let us now consider the counterterm for the momentum statistics. We will parametrize the time dependence of the 
counterterm as a power law in linear growth D as qD'^ and with S{k, a) = (fe, a) -|- QD'^{a)S^^^ (fe, a) we 

have then 

d^^\k,a) + (1 + j)cI qD'^{ k, a) 


fi{k, a) = S'{k, a) = fH 


(35) 
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FIG. 8 . Constraints on the speed of sound and its time dependence from the momentum-momentum auto power spectrum. 


In analogy to Eq. (30) we consider the estimator for the sound speed that first removes the one-loop corrections and 


maps the residual on the leading order counterterm 




- Pii - P 


Hoop 




= Cs,o>c^.0 9(2 + 7),c^,o(1 + 7) 


(36) 


The leading UV-sensitivity of SPT suggests that Cg scales as cr^, i.e., as (7 = 2). We saw above that this time 
dependence provided a decent description of the density power spectra. The constraints from the momentum-density 
and momentum-momentum statistic are shown in Figs. and respectively. As above, we again see that 7 = 2 
performs very well both on large scales, but the errors are considerable and Cg constraints for the momentum statistics 
show a strong scale dependence, as they did for the density. 

We will thus continue to consider the the 7 = 2 scale dependence and calculate the two-loop momentum coun¬ 
terterms. The counterterm for the density-momentum correlator can be obtained from the time derivative of the 
density-density counterterm or equivalently from taking the respective limits of Eq. (32) 


PSfj,,ctr — 




■3(2PS’ 




pgi—)-oo\ 
■ ^ 33-11 ) 




(37) 


We see in Fig. that this counterterm in combination with the finite part of the two loop calculation can indeed 
describe the scale dependence of the Cg estimator for k > 0.07 ZiMpc”^. The SPT inspired time dependence seems 
to work for this scale dependent part up to z = 1, but at z = 2 there are some deviations at fc > 0.2 /iMpc”^. We 
also overplot the effect of IR-resummation on the two-loop result, which will be described in more detail in the next 
section. 

To obtain the expression for the momentum-momentum correlator one would need to write the expressions for the 
two-loop counterterms at the field level and take the time derivatives. Here we rather follow the simple approach and 
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FIG. 9. Effect of the IR-resummation on the one- and two-loop power spectra. We show the ratio of the power spectra before 
and after IR-resummation with respect to the corresponding no-wiggle power spectrum in order to remove broadband effects. 
Below k — 0.2 /iMpc“'^ the bare two-loop calculation agrees with the IR-resummed one loop calculation at the percent level. 
After IR-resummation the wiggle part of the one- and two-loop calculation agree, which tells us that the IR-resummation 
captured the relevant terms in the explicit two loop calculation correctly. 


consider the limits of Eq. (331 to obtain 


— .P/^/i,ctr,lloop 


I0p«l^oo 




qpqi^oo 


P 


/j./i,ctr,21oop 


(38) 


The Cg constraints from this statistic are shown in Fig. While the agreement between the scale dependence of the 
two-loop calculation and the data is somewhat worse than for the density and density-momentum cross correlation, 
the model is still able to roughly capture the broadband scale dependence and the wiggles in the Cg estimator. 


B. IR resumation 

The EFT corrections discussed so far mainly address the broadband, i.e., short scale behaviour of the power spectra. 
Another set of corrections that are not fully captured by SPT are the long wavelength motions, which however do 
not affect the broadband behaviour due to the equivalence principle. The most prominent effect of the long modes is 
to damp the BAO oscillations. Lagrangian Perturbation Theory captures the effects of these motions better, since it 
keeps them resummed in an exponential and thus captures their effects to higher orders than the one to which the 
displacement field has been calculated. SPT only keeps motions to the order explicitly considered in the calculation. 
An effective way to combine the merits of both approaches is the so called infrared-resummation (IR-resummation) 
M: which calculates the broadband in SPT and corrects the result to account for the IR-motions. 

In this section we implement the simple IR-resummation described in |38j . This method multiplies the oscillatory part 
of the power spectrum (the wiggle part Pw) by an exponential damping but leaves the broadband part (the no-wiggle 
part Pnw) unaffected 




(l + Eg(fc)fc^) Pll,w + Plloop.w 


Pi 


lloop,nw 


(39) 


where 


S^(fc) = ^ / -p(<?fBAo)+2j2(gr-BAo)] , (40) 

with j„ being the n**' order spherical Bessel function. We only consider the smoothing due to motions arising from 
scales much larger than the scale under consideration, for definiteness we choose T,g{k) = T,{k/2). 

The explicit-two loop SPT calculation contains the resummed motions up to second order, such that the resumma- 
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FIG. 10. Ratio of the data to the various PT models at redshifts 2 = 0,0.5,1,2 from top left to bottom right. We show the 
linear theory calculation (green dot-dashed), the one-loop EFT (red solid) and the two-loop calculation (blue dashed). For the 
EFT calculation we show results both before (thin) and after IR-resummation (thick). The ratio is evaluated at the simulation 
data points and the two sigma errors on these data points are indicated by the gray band. 


tion needs to start from third order, i.e., E®. 


Pm (^1 -b E^(fc)fc2 + iE4(fc)fc4^ Pn.w + (l + ^iioop.w + -P 2100 P.V 

~b “b .Plloop,nw “b ^21oop,nw ■ 


(41) 


The effects of IR resummation are highlighted in Fig. Performing the IR resummation on the bare one-loop 
calculation leads to considerable changes to the power spectrum. Below k « 0.2 /iMpc”^, the not IR-resummed 
two-loop calculation performs almost as well as the IR-resummed one loop calculation. The IR-resummation of the 
two-loop calculation only matters at the percent level for k > 0.2 /iMpc”^. 

As we have seen above in Fig.|^ the two-loop calculation is tracking part of the BAO wiggles in the power spectrum 
residuals after the one-loop result has been removed. Let us now study its performance at higher wavenumbers and 
in the power spectrum itself. In Fig. we show the performance of the IR-resummed and not IR-resummed one- 
and two-loop EFT calculations with respect to the non-linear power spectrum extracted from the A^-body simulation. 
Let us first discuss the broadband performance. At redshift z = 0 the one loop calculation extends the range of 
validitj0of linear theory from k « 0.05 /iMpc”^ to fc « O.I ZiMpc”^. This is significantly less, than usually considered 
for the range of validity of the EFT at redshift z = 0 and arises from the fact that we have fixed the leading order 
counterterm in a way that is compatible with the largest available scales. We then use this parameter to calculate 
the two loop counterterm. This term, together with the finite part of the regularized two-loop calculation allows us 


® For the sake of definiteness we will commonly consider 1% deviations from the theory as the threshold for the range of validity. Many 
applications will require tighter errorbars on large scales to fix the amplitude. On smaller scales we will anyways suffer from baryonic 
effects and significant covariance, such that less restrictive requirements could be employed. 
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to extend the 1% agreement range to fc « 0.3 Here we should stop for a minute and reconsider the goal 

of this exercise. Usually one tries to fit the non-linear power spectrum as well as possible up to the highest possible 
wavenumber. But actually this should not be the goal of the fit with the deterministic part of the EFT, which we 
are computing here. The non-linear power spectrum is the sum of this deterministic part and the stochastic part. As 
we have pointed out in mi [32], this stochastic term amounts to a percent of the total power at fc = 0.25 /iMpc ^ 
(and about 3% at A: = 0.3 AiMpc”^). This means that the deterministic part of the power spectrum should deviate 
from the non-linear power spectrum by at least this much for k > 0.25 /iMpc”^. The deterministic EFT calculation 
(performed here) should asymptote to the perturbative/deterministic part of the power spectrum Ppx = Pni — J^stoch 
and not to the non-linear power spectrum itself. Thus, we have slightly underfitted the Cg parameter, overfitted the 
power. Once our EFT calculation is failing, it predicts more power than the non-linear power spectrum (downturn 
in Fig. 10). This failure would have happened even earlier and more steeply if the (positive) stochastic contribution 
had been subtracted from the non-linear power spectrum. A slight increase of Cg would bring our curve closer to the 
deterministic part of the power spectrum. Furthermore, looking at Fig. [^ , we see that both the one- and two-loop 
terms and their counterterms act to decrease power for a while before crossing zero and adding power. We might thus 
expect the higher loop calculation to reduce power in the same way at fc > 0.25 /iMpc”^ such that the prediction 
matches the deterministic part of the power spectrum for these wavenumbers. 

Having said this, let us now consider the effect of IR-resummation. As we have seen before, the two-loop results before 
and after IR resummation agree for fc < 0.2 fcMpc”^. For higher wavenumbers the IR-resummation indeed reduces 
the amplitude of the residual BAO wiggles even up to scales where the EFT broadband significantly deviates from 
the non-linear power. 


IV. OUTLOOK TO HIGHER ORDERS 

While the range over which the perturbative or deterministic part of the EFT can describe the full non-linear 
structures is certainly limited by stochastic terms, there might be some hope for further leverage in either precision 
at low wavenumbers or reach at intermediate wavenumbers from going to higher orders for the deterministic part. As 
shown in |7], a three-loop calculation is in principle feasible, but not very useful when taken at face value. As we have 
seen above, part of the higher loop calculation is degenerate with lower order counterterms. At the three-loop level 
this amounts to identify the terms that have already been taken care of by the terms that lead to either Puoop.ctr or 

.fc21oop,ctr- 

The three-loop power spectrum is given by 


^31oop — 2 Pi7 -|- 2P26 + 2P35_/ -|- 2P35_77 -|- Pm-I + ^44-77 , 


(42) 


with constituent diagrams given in Fig. El 
Let us discuss these terms separately: 

• P 17 : This is a propagator term. Its leading order contribution for low external wavenumbers scales as k^Pn 
and is thus degenerate with the leading order counterterm. To the extent that Cg has been fixed, we can 
safely subtract this term from the three-loop calculation. We expect such a contribution from the triple hard 
limit. Then there will be double hard limits, where one of the loops remains at the same order as the external 
momentum. These terms have been accounted for in Pis^ctr and thus need to be subtracted from the result. 
The part that is new, is the one where only one of the loops becomes large and this limit will become part of 
the three-loop counterterm. 

• P 2 e- The double hard limit of this term is degenerate with P 24 ,ctr and thus needs to be subtracted from the 
calculation. The single hard limit will lead to a new counterterm of the form fc^ x finite. 

• P 35 _ 7 : The nature of this term is very similar to P 24 , just extended by one loop of the stochastic kind. There are 
no degeneracies with lower order counterterms but it leads to a new two-loop counterterm of the form fc^ x finite. 

• P35-11 This term is a product of Pn, P13/P11 and P 15 /P 11 . The fc^Pn part of the P 15 contribution needs to 
be subtracted from the diagram. Finally, the counterterm is a combination of Cgfc^As and Cgfc^As.ctr 

• P44-1 is of the same nature as P33-7 and P22, be., the limits of high loop momenta lead to stochastic terms and 
the amplitudes of the divergencies are suppressed by 9^. 

• P44-11 The double hard limit with both “ear” diagrams large is of the form fc'^cr)^ x finite and thus leads to a 
new counterterm. 
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Pn P35-1 Pii-i 




FIG. 11. Diagrams for the three loop calculation. 


In summary, we can conjecture for general Z-loop diagrams that the diagram in which all the legs of the kernels 
are connected to another are stochastic terms and will not contribute to the leading order counterterm. All 

the other diagrams are dressings of diagrams that were encountered at lower orders with ’’ears”, i.e. diagrams that 
are given by the pair of momenta (< 7 , —q) in the kernels. When one of them is hard, we need a new counterterm at 
leading order in derivatives and when several of them are hard, this term will be degenerate with a counterterm that 
has already been fixed. Thus their contribution has to be subtracted from the finite part of the loop calculation. 


V. CONCLUSIONS 

Our study leads to the following conclusions: 

• Fitting the one-loop counterterm Cg at fc = 0.2 /iMpc”^ overestimates this coefficient by 50% for lower wavenum¬ 
bers and thus leads to too much suppression on large scales. While this regime is hard to extract from simulations, 
this might have important consequences for studies of primordial non-Gaussianity, where both the maximum 
wavenumber but also precision on large scales matter. 

• We find that the finite part of the two-loop calculation that is not degenerate with the leading order counterterm 
induces a 1% correction in the power spectrum at fc = 0.1 /iMpc”^ for our ACDM cosmology at redshift z = 0. 
The two-loop counterterm itself adds another 1% correction in the same direction at fc = 0.16 /iMpc”^. 

• We have not explicitly calculated the three-loop correction, but based on a related study in Lagrangian space 
m and its extension to density fields [35], we estimate that the stochastic term will lead to a 1% (3%) correction 
at fc « 0.25 /iMpc”^ (k « 0.3 /iMpc”^). Note that the perturbative EFT should describe the deterministic part 
of the power spectrum, which by definition has less power than the full non-linear power spectrum. Our EFT 
model fails by overpredicting the power compared to simulations so the failure of the model really happens at a 
slightly lower wavenumber. Also, there does not appear to be much room for the explicit three loop calculation 
to improve on our results without explicitly modelling the stochastic term. 

• We find that the two-loop EFT model with IR resummation can capture the scale dependence of the matter 
power spectrum up to fc « 0.3 ZiMpc”^ at the 1% level. Besides that, with the same EFT parameter and our 
assumption of a particular (D^) time dependence of the counterterms, we are also able to explain the scale 
dependence of the momentum power spectrum. 

• The fact that the EFT model can match the power spectrum does not necessarily mean that the EFT is the 
right model for the density field. The model so far does not include the stochastic term, that is partially given 
by the virialized motions within haloes. We do expect the one halo term to play a role at the percent level 
around k « 0.3 ft,Mpc~^ based on our recent study |32|. The one halo term is just one contribution to the 
stochastic term. 

• At the sub-percent level, many numerical effects can affect the agreement between theory and simulations. On 
the one hand, one can use the perturbative calculations on very large scales to check linear growth in simulations. 
On the other hand, precise comparisons of simulations and theory on large scales might require the theory to 
be evaluated for the same seeds that were used to initialize the simulations as we discuss in App. 
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• In summary, the picture at this point looks fairly consistent, but it is certainly too early to claim hnal success. 
Measurements of Cg without cosmic variance on large scales as discussed in Appendix ^ would certainly be the 
cleanest way to get the leading order EFT correction as well as its time dependence. Unfortunately there seem 
to be convergence issues in the A^-body simulation. Besides that there might be several other effects that might 
warrant further study. For example, the assumption of the EdS correspondence for the growth factor might 
impact our results at the percent level once the one-loop corrections become order unity corrections [39| . 
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Appendix A: Explicit Formulae 

In this Appendix we write down the explicit formulae for the one- and two-loop contributions to set the normalization 
for the terms used in the main text. The corresponding diagrams are given in Fig.[^ For the one-loop terms we have 


Pl3=mk) [ F3{qi,-qi,k)P{qi) , 

J qi 

[ \F2{qi,k-qi)\'^P{qi)P{\k-qi\) , 

J Oi 


(Al) 


P 22 =2 


where the recursion formulae for the SPT kernels can be found in the literature [ 3 ]. The kernels employed here are the 
ones derived for an Einstein de-Sitter (EdS) Universe and are based on separability of spatial and temporal structure 
of the theory. Note that all the kernels are symmetrized over their arguments. Their time dependence is governed by 
the scaling of the linear power spectrum with the growth factor D. Eor the two-loop terms we have 

Pi5=15P{k) f j F5{qi,-qi,q2,-q2,k)P{qi)P{q2) , 

Jqi Jq2 
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Note that we have not symmetrized the Pij terms for i ^ j, i.e., their contribution to the total equal time power 
spectrum will be 2 Pij. 


Appendix B: Cosmic variance and the propagator 

The propagator [ 4 ] measures the response of the density field to the initial conditions 

- ,')G(k) = (1^) , (Bl) 

It can be estimated from the cross power spectrum of the non-linear and the linear held and the auto power spectrum 
of the linear held as 


G{k) = . (B 2 ) 

Pii{k) 

The propagator is in principle the cleanest observable for the extraction of the speed of sound Cg of the Eulerian 
EFT, since the latter is nothing but a modihcation of the low-fc limit of the propagator. At one-loop level we have 
G = 1 -I- P13/P11 — c^k'^ and at two loop level G = 1 -I- P13/P11 — c^k'^ + P15/P11 + As.ctr/Ai- Here, we again use 
the regularized P13, from which all the k^P contributions have been removed. This is consistent with the previously 
introduced logic, where Cg absorbs the diherence after the one loop calculation has been subtracted from the data 
and effectively contains all the higher order and non-perturbative corrections as well as the mistake in the one-loop 
calculation. As argued before, this choice is somewhat arbitrary. This is no problem, as long as this choice is employed 
consistently. RPT [i] tries to resum the propagator in the high-fc limit, whereas here our goal will be to estimate the 
leading order corrections in the low-fc limit. 

We are interested in the deviations of the propagator from unity and define an alternative estimator for Cg 

Hni,i — Pii — P13 „ _ Pn\,i — Pii — P13 — P15 

PPii ’ “ PPii 
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FIG. 12. Left panel: Speed of sound extracted from propagator measurements in our simulations. We show measurements 
for two different choices of the Gadget PMGRID parameter PMGRID = Np and PMGRID = 2Np for the L simulation as 
well as for the M simulation with PMGRID = 2Np. The upper panel shows the k^P coefficient after the subtraction of the 
linear and one-loop corrections and the lower panel after subtraction of the regularized two-loop contribution (which are 
shown as solid and dashed lines in the upper panel). Right panel: Comparison of the estimator based on an analytical 
one-loop calculation and the one-loop calculation performed on the simulation IC grid. Once corrected for 3 x 10^ growth 
factor normalization the low-A: limit of the estimator in the simulations indeed asymptotes to a constant. 


The SPT contributions in these expression are IR sensitive, since the cancellation of IR modes in 2 Pi 3 -|- P 22 or 
its two-loop equivalent are not happening. When subtracting the one-loop SPT contribution from the propagator 
measured in simulations we address this issue by evaluating perturbation theory on the initial condition grid employed 
for the simulations using a technique similar to |40j . Thus, we are using the same IR modes that affect the non-linear 
dynamics in the simulation, thus directly addressing the IR sensitivity. The propagator based estimator projects out 
all the terms from the field that correlate with the linear field and thus provides an alternative and in some sense 
cleaner measurement of corrections that have the form of the leading order EFT counterterm, which could be masked 
by other contributions in the auto-power spectrum. 

We show the measurements of Cg based on the propagator in Fig. 


shows Cg before the two-loop contribution has been subtracted. 


12 


at one- and two-loop level. The upper panel 
Once corrected for 10“^^ level offsets in the linear 


growth factor]^ that are probably related to the timestepping in GADGET and that lead to an upturn or downturn 
of the data points in the plot, we see that the data asymptote to a constant on large scales and then decay on 
smaller scales. The shape of this decay is however captured by the scale dependence of .P 15 . Another remarkable 
observation is that the measured value of Cg depends very strongly on the PMGRID parameter in GADGET, leading 
to a shift of Acg « 0.8 /i“^Mpc^. We had seen a similar sensitivity already for the Lagrangian EFT coefficient of 
the displacement field in (31]. Based on this study, we are inclined to favor the results of the PMGRID = 2Np case. 
For another observation of this sensitivity in the power spectrum see |41j . Further evidence for the trustworthiness 
of this case comes from the fact that it agrees with the results from the higher resolution box, the M simulation. In 
a second step we now remove the scale dependence of .P 15 and see in the lower panel, that the estimated Cg is flat up 
to k K 0.2 /iMpc”^. There is clearly more need for convergence studies of the propagator and we certainly do not 
want to overinterpret a result that is so sensitive on numerical parameters of simulations. We conclude however, that 
there is evidence for a non-zero speed of sound correction after one loop SPT has been subtracted. The inferred value 
Cg = 1.15 /iMpc”^ roughly agrees with the value employed for the equal time correlators in the main text. Had we 
subtracted the explicit low-fc limit of P 15 (not P 15 , i.e. before regularization), the estimated Cg would have changed 
by Acg = —2.58 ZiMpc”^ and thus yielded a non-zero, negative Cg. The explicit corrections from three loops are even 
higher [7]. 

we also compare the grid based calculation for the Cg constraint from the auto power 


12 


In the right panel of Fig. 

spectrum with the analytical calculation. For the latter, we saw in the main text in Fig. [^that the data points for Cg at 
k < 0.05 /iMpc”^ are systematically low. This problem vanishes once the theory is calculated on the simulation grid. 
We see that the Cg estimator now asymptotes to a constant horizontal line on large scales, as we would expect it to 


We have independent evidence for such an error at this level from comparing power spectra for our fiducial parameter settings with a 
simulation with smaller timesteps. 
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FIG. 13. Low-fc measurements of from the final field and the propagator after the two-loop terms have been corrected 
for. The dashed line shows the effect of a relative error of the linear growth factor of 2 x 10“^ that both statistics have been 
corrected for. 
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TABLE I. Table of the two loop limits, references to the equations where they are discussed, the power of the cutoff dependence 
for a power law power spectrum P{k) oc fc"' with general power law slope n and for n = —3/2. For the single hard limit the 
slope gives the power of the hard integral ignoring the remaining finite integral, while for the double hard integrals we consider 
both momenta in the hard integrals to be of the same order. The choice n = —3/2 is motivated by the slope of our ACDM 
power spectrum at fe ~ 0.1 hMpc“^. 


based on the scale dependence of the two-loop corrections shown in Fig.[^ Again, the value of this asymptotic constant 
depends strongly on the PMGRID parameter choice, now leading to a Acg = 0.3 /i“^Mpc^ difference between the 
two cases. Note however that they agree at higher wavenumbers. Thus, to the extend that our ansatz is trustworthy, 
a model what matches at these scales would prefer the PMGRID = 2Np case at lower wavenumbers. 

There is also a slight disagreement between the propagator and power spectrum estimates for the favored PMGRID= 


2Np case. The power spectrum method of this case would indicate a cf = 1.05 h ^Mpc . In Fig. 


13 


we show both 


the propagator and the power spectrum estimator after the finite two loop terms have been subtracted out. Except 
for a Acg « 0.1 h~^Mpc^ offset both estimators are flat and consistent up to /c « 0.15 hMpc”^, where higher order 
terms, for instance the two loop counterterms, start to matter. 


Appendix C: Limits of the two loop terms 

In the main text, we have concentrated our discussion on the terms that we consider relevant for the leading UV 
sensitivity and the corresponding counterterms. Let us, for the sake of completeness, discuss the remaining hard limits 
in this appendix. An overview of all the single- and double-hard limits of the two loop calculation is given in Tab. [l] 
In this table we also give the power of the cutoff dependence of the remaining integrals if the initial power spectrum 
is of power law form P{k) oc fc". We evaluate the cutoff dependence for n = —3/2, the slope of our power spectrum 
at fc = 0.1 /iMpc”^. For the single-hard limits we immediately see that the terms that we found to dominate the 
shell behaviour have the most shallow decay in the UV, and are thus the most sensitive to the change of the power 
spectrum at high wavenumbers. For the double hard limits, the limit of P 15 is still growing for n = —3/2 but turns 
around at for n = —2 at fc « 0.3, so it will still converge based on the high-fc slope of our initial power spectrum. Yet, 
it is immediately clear why this integral should be absorbed into the counterterm. The subleading k^P UV-sensitivity 
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of Pi 5 (not mentioned in the table but below in Eq. C8) scales as 2n + 2, i.e., as for n = —3/2 and should thus be 
the next term considered as a counterterm, after the single hard limits. It will change the coefhcient of the 
counterterm. 

First, we consider the limit of P33-1 for <71 —>■ 00 while <72 remains finite 

32 ^ r 5565rf - 11465r| + 26409r| + 22285r| - SlOrz 


J J Qo 'I 


14288400rf 


15(ri-l)^(7ri + 2) (53ri + 17) / l + rz 

^ V 1 - "^2 


14288400rf 


(Cl) 


The amplitude is given by / Pii{q)/q'^, which is quickly convergent in the UV and does thus not contribute to a 
significant UV sensitivity. For the double-hard limit of the same integral we have 


P33-I^°° f [ K33-l{QljQ2)Pll{qi)Pll{Q2)P{\Ql + (l2\) , 
Jqi Jq2 


(C2) 


where K 33-1 can be parametrized in terms of the magnitude of the momenta and their cosine /ii 2 = <71 • < 72/9192 as 


K33-1 = 


(Mi 2 — 1) (9192 ( 392 / 1^2 T 6347x12) -f (3927X12 -f 121) Qi -|- (3927X12 -f 121) (7I) 
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The term thus scales as in the fc —>■ 0 limit, with a quickly convergent amplitude J J P/1/9^. 
For the double-hard limit of P 24 we have 


(C3) 


091.92 ->-00 

^24 


--k^ f f ^ 24 ( 91 , 92 )^ 11 ( 91 )^ 1 ( 92 )^ 11 ( 191 + 921 ), 

Ja^ JQ'> 


(C4) 


where K 24 is given by 

^24 = 
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The term thus scales as in the fc —>■ 0 limit, with a quickly convergent amplitude. 

In the main text we considered only one of the single-hard limits of P 24 , the one where the closed loop or “ear” 
diagram becomes large. The other limit, where the loop containing two power spectra becomes large is given by 
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with a kernel 
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As for the amplitude is given by / Pli{q)/q'^. 
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As we had hinted at in the main text, if both loop momenta in P 15 go to infinity, we end up with an analytic 
function in k 

P/5'’9=^~ =-fc2Pn(fc) / / 7 ^ 15 . 2 ( 91 ,92)Pii(9i)7^ii(92)+fc^Pii(fc) / / 7 ^ 15 , 4 ( 91 ,92)Pii(9i)7^ii(92) - (C 8 ) 

Jqi Jq2 Jqi Jq2 

We rederive the double hard limit of P 15 and note that the expression agrees with what was found for the k^Pn 
part in 
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Note that the limit of the above expression, where one of the momenta is smaller than the other one <72 ^ 9i, turns 


out to be the same as the limit of Eq. (17), in which we first take qi ^ k,q 2 and only afterwards q 2 ^ k. We also 
obtain the expression for the k'^Pn part to estimate the size of the subleading UV-sensitivity 
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